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ABSTRACT 

The recent detection of a 3.5 keV X-ray line from the centres of galaxies and clusters by 
Bulbul et al. (2014a) and Boyarsky et al. (2014a) has been interpreted as emission from the 
decay of 7 keV sterile neutrinos which could make up the (warm) dark matter (WDM). As 
part of the COpernicus COmplexio (COCO) programme, we investigate the properties of dark 
matter haloes formed in a high-resolution cosmological A^-body simulation from initial con¬ 
ditions similar to those expected in a universe in which the dark matter consists of 7 keV 
sterile neutrinos. This simulation and its cold dark matter (CDM) counterpart have ^ 13.4bn 
particles, each of mass ^ 10^ h~^ Mq, providing detailed information about halo structure 
and evolution down to dwarf galaxy mass scales. Non-linear structure formation on small 
scales (M 200 ^ 2 X 10^ h~^ Mq) begins slightly later in COCO-WARM than in COCO-COLD. 
The halo mass function at the present day in the WDM model begins to drop below its CDM 
counterpart at a mass ^ 2 x 10^ h~^ Mq and declines very rapidly towards lower masses 
so that there are five times fewer haloes of mass M 200 = 10^ h~^ Mq in COCO-WARM than 
in COCO-COLD. Halo concentrations on dwarf galaxy scales are correspondingly smaller in 
COCO-WARM, and we provide a simple functional form that describes its evolution with red- 
shift. The shapes of haloes are similar in the two cases, but the smallest haloes in COCO-WARM 
rotate slightly more slowly than their CDM counterparts. 
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1 INTRODUCTION 

The identity of dark matter, the dominant matter component of the 
Universe, has long been a subject of great interest in cosmology. 
In the last three decades, the model of non-relativistic dark mat¬ 
ter consisting of heavy weakly-interacting particles with negligi¬ 
ble thermal velocities at early times, the Cold Dark Matter (CDM) 
model, has become the cornerstone of the standard cosmological 
paradigm. The standard model with dark energy in the form of a 
cosmological constant, A (ACDM, henceforth just CDM) has been 
very successful in predicting and matching observational data on a 
wide range of scales, from the temperature fluctuations in the Cos¬ 
mic Microwave Background (Planck Collaboration et al. 2014) to 
the statistics of galaxy clustering (Colless et al. 2001; Zehavi et al. 
2002; Hawkins et al. 2003; Tegmark et al. 2004; Cole et al. 2005; 
Eisenstein et al. 2005; for a comprehensive review on the subject, 
see Frenk & White 2012). 

* E-mail: sownak.bose@durham.ac.uk 


With the advent of the LHC it was hoped that one of the best- 
motivated CDM candidates, the lightest supersymmetric particle 
(the neutralino) would be found. The lack of evidence for super- 
symmetry at the LHC and the absence of a convincing direct or 
indirect signal for CDM (but see Hooper & Goodenough 2011) has 
encouraged the exploration of viable alternatives. One of the most 
promising alternatives is the sterile neutrino (Dodelson & Widrow 
1994; Asaka & Shaposhnikov 2005), which behaves as warm dark 
matter (WDM) due to the particles’ non-negligible thermal veloci¬ 
ties at early times. Being collisionless, this leads to free streaming 
and the damping of perturbations in the density field, creating a 
cutoff in the matter power spectrum on the scale of dwarf galaxies. 

A simple extension of the Standard Model of particle physics, 
called the neutrino Minimal Standard Model (z/MSM, Boyarsky 
et al. 2009), consists of three right-handed sterile neutrinos in 
which, for a specific choice of parameters, one of the sterile neutri¬ 
nos behaves as a dark matter particle and the model explains neu¬ 
trino flavour oscillations. Each one of this triplet of particles has its 
mass below the electroweak scale; one in the keV scale (denoted by 
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Ml), and two in the GeV scale (denoted by M 2 and M 3 ). The for¬ 
mer behaves as a relativistic particle at the time of neutrino decou¬ 
pling and acts as WDM, and is then redshifted to non-relativistic 
energies during the radiation-dominated era. Unlike a thermal relic, 
the cutoff in the power spectrum introduced by a sterile neutrino 
of a fixed mass depends on a second parameter, the lepton asym¬ 
metry. As we explain later in the following section, it is possible to 
approximate the sterile neutrino power spectrum with a WDM ther¬ 
mal relic equivalent, particularly for very low and very high values 
of the lepton asymmetry. 

The unidentified 3.53 keV X-ray line originally detected in the 
spectrum of a stack of galaxy clusters (Bulbul et al. 2014a) and in 
the spectra of M31 and the Perseus cluster (Boyarsky et al. 2014a) 
could be a decay signal of sterile neutrino dark matter, with a par¬ 
ticle mass of 7 keV. More recently, Boyarsky et al. (2014b) have 
also identified a similar line in the centre of the Milky Way. There 
are, however, several groups that have questioned the plausibility 
of this detection. For example, Riemer-Sorensen (2014) failed to 
find a signal in Chandra observations of the Milky Way. Of course, 
the Galactic centre is heavily contaminated by X-rays, which intro¬ 
duces uncertainties, a point made by Boyarsky et al. (2014b). 

Systematic effects can result from the atomic data used in 
modelling the plasma, as argued by Jeltema & Profumo (2014), 
who found no excess when re-analysing the Boyarsky et al. (2014a) 
data and claimed that any signal at 3.5 keV could be explained by 
known Potassium (K XVIII) and Chlorine (Cl XVII) lines. Bulbul 
et al. (2014b) put this latter result down to the use of “incorrect 
atomic data and inconsistent spectroscopic modelling” by Jeltema 
& Profumo (2014). A further non-detection was then reported in 
the stacked spectra of galaxies from Chandra and XMM-Newton 
(Anderson et al. 2014), while most recently, Malyshev et al. (2014) 
analysed the spectra of stacked dwarf galaxies from XMM-Newton 
and claimed to rule out the Andromeda signal detected by Bulbul 
et al. (2014a) at the 4.6cr confidence level. This has spurred other 
groups (see, for example, Conlon & Day 2014) to associate the 3.53 
keV signals to the conversion of a sterile neutrino into an axion, 
and its subsequent decay into photons. Such a mechanism requires 
a magnetic field, the presence and strength of which can vary from 
galaxy to galaxy, a scenario that could explain why this line is only 
seen in some objects. 

Clearly, whether or not the 3.53 keV line really does exist re¬ 
mains an open question. It is, therefore, important to investigate 
the predictions for the formation of cosmic structures in a model 
in which the dark matter consists of particles that could decay pro¬ 
ducing such a line. Constraints on such models can be set from the 
observed clustering of the Lyman-a forest at high redshift whose 
small-scale structure would be erased if the dark matter were warm. 
On these grounds, Viel et al. (2013) recently set a (current) lower 
limit of 3.3 keV for the mass of a dominant thermal warm dark 
matter particle. 

Coincidentally, the power spectrum of a 3.3 keV thermal 
warm dark matter particle is well approximated by that of a 7 keV 
sterile neutrino for a lepton asymmetry of Lq = 8 . 66 . This corre¬ 
sponds to the smallest allowed value of the power spectrum cutoff 
length (i.e. to the “coldest” power spectrum possible) for a sterile 
neutrino of mass 7 keV. This is the model that we will explore in 
this work. Ruling out this model from astronomical data on small 
scales would rule out the entire family of 7 keV sterile neutrino 
candidates. To investigate the model we use high resolution N- 
body simulations whose results we compare with those of CDM 
simulations with the same phases in the initial conditions. We are 
interested exclusively in characterising the properties of dark mat¬ 


ter haloes of mass in the region of the power spectrum cutoff and, 
in this study, we ignore the effects of baryons. Such effects must be 
taken into account when comparing model predictions with obser¬ 
vations. In the case of CDM, relevant baryon effects on the small 
scales of interest here have recently been quantified by Sawala et al. 
2013, 2014, 2015; Schaller et al. 2015) 

The layout of this paper is as follows. In Section 2, we intro¬ 
duce the simulations used in this work, the modelling of the WDM 
component, and describe how we tackle the issue of spurious halo 
formation in our simulations. In Section 3 we present our main re¬ 
sults from the comparison of WDM and CDM from our simula¬ 
tions, in terms of both the large-scale distribution of matter, and 
the internal structure of haloes. Finally, in Section 4, we summarise 
our findings and look into some future work that will be carried out 
with the same set of simulations. 


2 THE SIMULATIONS 

In this section, we provide an overview of the initial conditions and 
modelling of the WDM component in our simulations. 

2.1 The simulation set-up 

The A^-body simulations presented in this paper are part of the 
COpernicus COmplexio (COCO) simulation programme (Hellwing 
et al. 2015) being carried out by the Virgo Consortium. This is a set 
of cosmological “zoom-in” simulations (Katz & White 1993; Frenk 
et al. 1996), as was done in the GIMIC simulations (Crain et al. 
2009). The parent simulation, called the COpernicus complexio 
LOw Resolution (or COLOR) simulation, followed the evolution of 
4.25 billion particles in a periodic box of size 70.4 Mpc. We 
extracted a roughly spherical region of radius ~ 18h~^ Mpc, and 
centred on the location (42.2, 51.2,8.8) h~^ Mpc in the COLOR 
volume. Both COLOR and COCO assume cosmological parame¬ 
ters derived from the seven-year Wilkinson Microwave Anisotropy 
Probe (WMAP 7) data (Komatsu et al. 201 1), with the parameters: 
Qm = 0.272, Qa = 0.728, h = 0.704, Us = 0.967 and erg = 
0.81. Here, 0.^rn,A} represents the present-day fractional contribu¬ 
tion of matter and the cosmological constant respectively, in units 
of the critical density pc = SHq/SuG, h = iTo/100 km/s/Mpc 
is the dimensionless Hubble parameter, Us is the spectral index of 
the primordial power spectrum, and ag is the linear rms density 
fiuctuation in a sphere of radius Sh~^ Mpc at z = 0. 

Dark matter particles with three different masses are used in 
regions simulated at different resolutions within the parent sim¬ 
ulation volume. Initially, the high-resolution region has a shape 
similar to an amoeba which approximates a sphere of radius ~ 
17.4h~^ Mpc at the present time. It contains 12.9 billion particles 
of mass 1.135 x 10^ h~^ Mq. The volume surrounding this region 
contains the medium- (3.07 x 10^ h~^ Mq) and low-resolution 
(1.96 X 10^h~^ Mq) particles. We have taken care to minimise 
contamination of the high-resolution region by lower mass particles 
and all the haloes discussed in this study are entirely made up of the 
high-resolution particles. The gravitational softening was kept fixed 
at e ~ 230 h~^ pc for the high-resolution particles, increasing by a 
factor of 10 each time for the medium- and low-resolution particles. 

The simulation ran from = 127 to = 0 using the GAD¬ 
GET 3 code, which is an updated version of the publicly available 
GADGET2 code (Springel et al. 2001a; Springel 2005). Phase in¬ 
formation for the creation of the initial conditions for both COCO¬ 
WARM and COCO-COLD was obtained from the public Gaussian 
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white noise field PANPHASIA (Jenkins 2013), and perturbations 
thereafter were calculated using the second-order Lagrangian Per¬ 
turbation Theory (2LPT) algorithm presented in Jenkins (2010). 
The details of the simulation, along with the PANPHASIA phase de¬ 
scriptor, are summarised in Table 1 . 

The distinctive feature of WDM particles are non-negligible 
thermal velocities at early times, which result in free streaming 
that washes out perturbations in the matter distribution below the 
free streaming scale (Bond & Szalay 1983; Schneider et al. 2012; 
Benson et al. 2013). As a result, we expect the abundance, distri¬ 
bution and internal structure of WDM haloes to be different from 
those of CDM haloes. Indeed, thermal velocities introduce a limit 
to the fine-grained phase space density in dark matter haloes, cre¬ 
ating cores in the density profile (Maccio et al. 2012; Shao et al. 
2013). However, as shown in these papers, the cores produced by 
realistic thermal relics are only a few parsecs in size, and thus not 
astrophysically relevant. In our simulations we can neglect these 
thermal velocities, which at z = 0 are of the order of a few tens 
of metres per second (Lovell et al. 2012) so, over the course of the 
simulation, which starts diiz — 127, the particles would travel only 
a few kiloparsecs, comparable to the mean interparticle spacing of 
the high-resolution particles. 

The WDM power spectrum of density fluctuations is often 
modelled by the transfer function, T(k), relative to the CDM case: 

Twdm(A:) = T^(A:)Tcdm(A^) • (1) 


We approximate T(k) using the fitting formula provided by Bode 
et al. (2001): 

T(fc)= + , (2) 


where a and o are constants. As computed by Viel et al. (2005), 
fork < 5 h~^ Mpc, the value u = 1.12 provides the best-fitting 
transfer function. The value of a is dependent on the mass of the 
WDM particle (Viel et al. 2005): 


r mwDM 1 

-1.11 

f^WDM 

0.11 

h 

1 keV J 


0.25 


0.7 


h ^ Mpc, 


(3) 


and determines the scale of the cutoff due to free streaming in the 
WDM power spectrum relative to CDM. It should be noted that this 
transfer function is a fit to the full thermal relic power spectrum, 
obtained by solving the Boltzmann equation (Lovell et al. in prep.). 

As we can see in Eq. 3, the “warmer” the dark matter particle 
(i.e., the lower its rest mass is), the larger the scale at which the 
cutoff in the power spectrum occurs. One way to define the char¬ 
acteristic scale in the power spectrum is through the “half-mode” 
wavenumber, /chm, where the transfer function in Eq. 2 drops by a 
factor of two: 


fchm = i (2"/" - 1) 


l/2iy 


(4) 


The associated “half-mode mass”, Mhm, is the mean density en¬ 
closed within this half-mode: 


7\ /f - f 


(5) 


Eor the 3.3 keV model, this occurs at around 
Mhm ~ 2 X 10^ h~^ Mq (Colin et al. 2008; Angulo 
et al. 2013; Viel et al. 2013). We will show later that differences in 
the formation time of haloes in WDM and CDM begin to appear 
below ~ 2 X 10^ h~^ Mq, approximately an order of magnitude 
above the half-mode mass scale. 


The form of the power spectrum for sterile neutrino dark mat¬ 
ter is determined by two parameters - the mass of the particle, 171 ^^, 
and the lepton asymmetry, Lq : 

Lg = 106 , (6) 

where rijy^ is the number density in electron neutrinos, the 
number density in electron anti-neutrinos, and s is the entropy den¬ 
sity of the Universe (Laine & Shaposhnikov 2008). The scales at 
which the power spectrum is suppressed for sterile neutrinos vary 
non-monotonically as a function of Lq. If Lq is very small (<^ 1) 
the power spectrum exhibits a similar abrupt cutoff to that of a ther¬ 
mal relic. As Lq is increased, the cutoff becomes gentler and khm 
shifts to larger values. At some value of Lq (typically between 8 
and 25 depending on the sterile neutrino mass), khm reaches a max¬ 
imum; for still higher Lq, khm retreats to lower k and returns to its 
original shape and position (Shi & Fuller 1999; Abazajian 2014, 
Lovell et al. in prep.). 

The power spectra used in the COCO simulations are shown as 
thick lines in Fig. 1 : CDM in black, 3.3 keV WDM in red and 7 keV 
sterile neutrinos with Lq = 8.66 in blue. All three power spectra 
agree on large scales. On small scales, the two warm dark matter 
models differ from CDM. khm for the sterile neutrino case occurs at 
a very similar scale, and the cutoff has a similar shape to that for the 
thermal relic case. On smaller scales still, the sterile neutrino power 
spectrum has more power than its thermal counterpart, but the dif¬ 
ferences only become significant on scales where the amplitude is, 
at most, a few percent of the peak amplitude. These differences are 
negligible and can be safely ignored in our simulations. The thin 
lines in the figure correspond to 7 keV sterile neutrino power spec¬ 
tra for different values of the lepton asymmetry, Lq. The Lq = 8.66 
model that we have simulated corresponds to the “coldest” possible 
7 keV sterile neutrino. 


2.2 Halo identification and matching 

Haloes were identified in our simulations using the friend-of-friend 
(EOF) algorithm (Davis et al. 1985) with a linking length of 0.2 
times the mean interparticle separation, and a minimum of 20 
particles. Gravitationally-bound substructures within these groups 
were then identified using the SUBFIND algorithm (Springel et al. 
2001 b), although in this paper, we will be mostly concerned with 
the properties of the WDM EOF groups. We determine the halo 
centre using the “shrinking sphere” method of Power et al. (2003). 
In short, we recursively compute the centre of mass of all parti¬ 
cles within a shrinking sphere, until a convergence criterion is met. 
In each iteration, the radius of the sphere is reduced by 5%, and 
stopped when only 1000 particles or 1% of the particles of the ini¬ 
tial sphere (whichever is smaller) are left. 

Comparing halo statistics between sets of simulations requires 
consistent definitions for the various properties of the haloes. In 
this work, we make use of two definitions of mass: Mfof, which 
is the mass of all particles identified by the algorithm as belonging 
to the EOF group, and M200, which is the mass contained within 
a sphere of radius 7200 (centred on the “shrinking sphere” centre 
defined above), within which the average density is 200 times the 
critical density of the Universe (pc) at the specified redshift. An¬ 
other common radius used to define a halo edge is the virial ra¬ 
dius, Tvir, within which the density of the halo p(< rvir) = Ap^ 
where A ~ 178U^^^ (motivated by the spherical collapse model. 
Eke et al. 1996. Note that this definition is consistent with the 
virial overdensity relation in Bryan & Norman 1998). Table 2 sum- 
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Parameter 

COLOR (Parent volume) 

coco (This paper) 

Box Size 

70.4 Mpc 


mwDM 

3.3 keV 

- 

Np 

4,251,528,000 

13,384,245,248 

Uhr 

70.4^ h-^Mpc^ 

2.2 X 10^ h-^Mpc® 

'^p,hr 

6.196 X 10® h-^ Mq 

1.135 X 10® h-^ Mq 

A^p,hr 

4,251,528,000 

12,876,807,168 

^hr 

1 

0 

230 pc 

h 

0.704 

- 

Um 

0.272 

- 

Ua 

0.728 

- 

0-8 

0.81 

- 

Phase descriptor 

[Panphl,L16,(31250,23438,39063),S12,CH582187950,COLOR] 



Table 1. Cosmological parameters used in the coco simulations, and its parent volume, color. Here, mwDM is the mass of the thermal relic warm dark 
matter particle, is the total number of particles (of all types) used in the simulation, 1/hr is the approximate volume of the high-resolution region at 2 ; = 0, 
^p,hr is the mass of an individual high-resolution dark matter particle, A^p^hr is the total number of particles of this species, whereas Chr is the softening 
length applied to them. The cosmological parameters /i, and ag are as described in the text. The phases for the parent color simulation can be 

generated using the panphasia phase descriptor provided in the last row. The blank fields in the coco column mean that the parameter assumes the same 
value as in the parent simulation, COLOR. 



log[A;/(AMpc ^)] 


Figure 1. The (dimensionless) matter power spectrum for: a thermal 3.3 
keV WDM (red), a sterile neutrino of mass miy^ = 7 keV and lep¬ 
ton asymmetry Lq = 8.66 (blue) and CDM (black). Both the WDM 
and sterile neutrino power spectra have significantly suppressed power at 
small scales, with the deviation from CDM case at almost identical scales: 
log (A:) > 1.0 A Mpc“^. Also shown as thin coloured lines are power spec¬ 
tra for 7 keV sterile neutrinos with different values of Lq , as indicated in 
the legend. 


Simulation Wfof(^ = 0) Arsubs(^ = 0) 


COLOR-COLD 

3,961,192 

4, 770, 041 

COLOR-WARM 

2,609,122 

3,082,275 

COCO-COLD 

8,896,811 

10,502,187 

COCO-WARM 

2,548,743 

2,830,514 


Table 2. Number of groups and subhaloes identified by the FOF algorithm 
and SUBFIND in COLOR and COCO at 2 ; = 0. 

marises the total number of groups and self-bound substructures 
identified at 2: = 0 in our simulations. 

Since both COCO-WARM and its COLD counterpart were sim¬ 
ulated using the same initial phases, we are able to match many 
objects between the two simulations. This allows us to make also 
object-by-object comparisons in addition to comparing just statisti¬ 
cal distributions of halo properties. In order to correctly match the 
haloes we do the following: first, we take the 50 most-bound parti¬ 
cles from a COCO-WARM halo, and look for the COCO-COLD halo 
in which there are at least 25 (50%) of these particles. We then 
confirm the match by repeating the same process, this time start¬ 
ing with the COCO-COLD haloes, in decreasing order of mass. This 
results in a bijective match between haloes in the two simulations. 
Using this method, we are able to match around 97% of haloes with 
M 200 > 10^ h ^ Mq- 

2.3 Spurious haloes and their removal 


Number counts of haloes and subhaloes are fundamental statistics 
of the halo population, so the correct identification of haloes is of 
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Figure 2. Redshift evolution of the projected dark matter density in COCO-COLD (left) and the 3.3 keV coco-warm Universe (right). From top to bottom, 
the top three panels show snapshots at 2 ; = 10, 2 : = 6, 2 ; = 1 of the projected mass density in cubes of side 2 h~^ Mpc, centred on the most massive group 
at 2 ; = 0. The bottom panels show zooms of a 5 x 10^° h~^ Mq halo at 2 ; = 0 in a cube of side 150 h~^ kpc. The emergence of small haloes at early 
times is apparent in the CDM case, when the WDM distribution is much smoother. The formation of large haloes occurs at roughly the same time in the two 
simulations and the subsequent growth of these haloes is similar in the two cases. In the zoom shown in the bottom panel, the lack of substructure in the WDM 
case compared to its CDM counterpart is stark. 
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primary importance. It has been known for some time (Wang & 
White 2007; Angulo et al. 2013; Lovell et al. 2014) that in simu¬ 
lations in which the initial power spectrum has a resolved cutoff, 
as is the case for COCO-WARM, small-scale structure is seeded in 
part by the discreteness of the particle set. In other words, a sub¬ 
structure finder will identify density peaks that have arisen not as 
a result of gravitational instabilities from a cosmological perturba¬ 
tion. These artificial fragments can often by identified “by eye” as 
they tend to be regularly spaced along filaments of the mass dis¬ 
tribution. They produce a power-law-like upturn at small masses 
in the WDM mass function. Since this is just a numerical (and 
resolution-dependent) artefact of our WDM simulations, care must 
be taken to identify these spurious haloes and, if appropriate, re¬ 
move them from the halo catalogue. While it is, in principle, pos¬ 
sible to eliminate these structures by increasing the resolution of 
the simulation, this is computationally prohibitive: Wang & White 
(2007) have shown that the mass at which spurious structures dom¬ 
inate the mass function scales with the number of particles in the 
simulation, asMocA^”^/^. 

Lovell et al. (2014) developed an algorithm to identify spuri¬ 
ous clumps in WDM simulations. A large number of them can be 
removed by performing a mass cut below a resolution-dependent 
scale, as suggested by Wang & White (2007): 

Miim = 10.Ip d , (7) 

where d is the mean interparticle separation and /cpeak is the spatial 
frequency at which the dimensionless power spectrum, A^(/c), has 
its maximum. Applying this condition on its own would also re¬ 
move some genuine haloes that form below this scale. Lovell et al. 
(2014) refined this criterion by also making a cut on the basis of the 
shapes of the initial Lagrangian regions from which WDM haloes 
form. They find that the spurious candidates tend to have much 
more flattened configurations in their (unperturbed) initial positions 
than genuine haloes, as judged from a CDM simulation. Defining 
the sphericity, s, of haloes as the axis ratio, c/a, of the minor to 
major axes in the diagonalised moment of inertia tensor of the ini¬ 
tial particle load, the sphericity cut is made such that 99% of the 
CDM haloes at that redshift lie above the threshold. 

Following exactly the methodology of Lovell et al. (2014), 
we clean the COCO-WARM catalogue as follows: (1) remove all 
(sub)haloes with Shaif-max < 0.165 irrespective of mass; (2) 
for those that pass (1), remove (sub)haloes with Mmax < O.fiMum. 
Here, Mmax is the maximum mass attained by a (sub)halo dur¬ 
ing its evolution, and Shaif-max is the sphericity (= c/a) of the 
(sub)halo at the half-maximum mass snapshot. This is chosen so 
as to identify a (sub)halo at a time well before it falls into a larger 
host, when its particles are subject to tidal stripping. The factor of 
0.5 in condition (2) is calibrated by matching between resolutions 
in the AQUARIUS simulations (see Lovell et al. 2014 for details). 
Having done so, we find that over 91% of the (FOF) haloes formed 
in COCO-WARM are in fact spurious, and are rejected from the halo 
catalogue when computing properties like mass functions. The ele¬ 
ments of this section are summarised in Fig. 3. 


^ The criterion Shaif-max < 0.165 is appropriate for haloes identified at 
2 ; = 0; for higher redshifts, one needs to determine the 1% sphericity cut at 
that redshift. 



Figure 3. Number density of haloes in the sphericity vs. maximum mass 
space in COCO-WARM at 2 ; = 0. The dashed black lines show the cuts 
on sphericity and mass that we use to clean the halo catalogue. Rejected 
(spurious) candidates are those that fail the cuts in the manner described in 
the text. 

3 RESULTS 

In both cold and warm models, dark matter haloes assemble in a 
hierarchical way, acquiring mass by merging with other haloes and 
by smoothly accreting ambient mass {e.g. Press & Schechter 1974; 
Frenk et al. 1985; Lacey & Cole 1993; Wechsler et al. 2002). In 
this section, we focus on global halo properties such as formation 
times, abundance and internal structure. We make a direct compar¬ 
ison between our cold and warm dark matter models. On scales 
much larger than the WDM suppression scale in the initial power 
spectrum, we expect the properties of haloes to be very similar in 
the two cases, but differences should become increasingly impor¬ 
tant at ~ 2 X 10^ h~^ M© and below. 


3.1 Redshift of formation 

The absence of primordial perturbations below the cutoff scale in 
the WDM power spectrum induces differences in the formation 
times of the smallest haloes. We can visualise these differences 
directly by examining the images displayed in Fig. 2. At early 
times, the projected density in COCO-WARM (right panels) is vis¬ 
ibly smoother than the equivalent projection in COCO-COLD (left 
panels), which has a “grainier” appearance owing to the very large 
number of haloes below ~ 10^ h~^ M© that form in this case, 
well before the first objects have collapsed in COCO-WARM. Thus, 
the onset of the structure formation process in this simulation is 
delayed relative to its CDM counterpart. 

In order to quantify the different halo formation epochs in 
COCO-WARM and COCO-COLD, we trace the evolution of each FOF 
group through its merger tree, and define the redshift of formation 
as the first time when the mass of the most massive progenitor ex- 
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ceeds half the final FOF mass: M (zform) = M ( 2 : = 0) /2 (e.g. 
Harker et al. 2006; Neto et al. 2007). Other definitions of halo for¬ 
mation time also exist in the literature (e.g. Navarro et al. 1996, 
1997), which should be borne in mind when making comparisons. 

The result, for all haloes in COCO-WARM (including spu¬ 
rious objects) and COCO-COLD is shown in Fig. 4. The forma¬ 
tion redshifts of haloes of mass M 200 ^ 2 x 10® h~^ Mq, are 
very similar in COCO-WARM and COCO-COLD, as expected. The 
difference between the two begins to manifest below a mass of 
M 200 ~ 2 X 10® h~^ Mq, an order of magnitude above the half¬ 
mode mass scale for a 3.3 keV WDM particle (c.f. Section 2.1). For 
these smaller haloes, 2;form is lower for WDM than CDM. The sud¬ 
den upturn in the WDM Zform for M 200 < 10^ Mq (shown in 
the open red circles) is a signature of the spurious haloes described 
in Section 2.3. From here on, we will exclude these spurious haloes 
and only show results from the cleaned COCO-WARM sample. The 
difference in formation times is a subject we will revisit when com¬ 
paring the concentration-mass relations of WDM and CDM in Sec¬ 
tion 3.4. Note that in this figure, we include all haloes, and not 
necessarily matched between CDM and WDM, which is why the 
medians at the largest mass bins are not exactly identical. 

We find that the delay in the formation time of CO CO-WARM 
haloes of a given mass, relative to COCO-COLD, is well described 
by the fitting function: 


WDM 

'^form 

^CDM 

^form 


^1 + a 


Mhm 

M 2 OO 


-b 


( 8 ) 


where Mhm is the half-mode mass introduced in Section 2.1, a = 
1.23 and b = 0.56. This fit is shown as the thin red line in Fig. 4. 


3.2 Differential halo mass functions 

Counting the number of dark matter haloes as a function of their 
mass is one of the simplest and most important population statistics 
that one can use to distinguish between WDM and CDM models, 
since fewer haloes will form in the former close to the half-mode 
mass. 

In Fig. 5, we show the build-up of the halo population as a 
function of redshift in COCO-COLD (solid lines) and COCO-WARM 
(dashed lines). The shaded regions and error bars represent the 
Poisson uncertainty in both cases. Spurious haloes have been omit¬ 
ted from the WDM differential halo mass function (dHMF) at each 
redshift, using the methodology outlined in Section 2.3. The edge 
of the grey region marks the nominal resolution limit of our simula¬ 
tion which corresponds to a halo with at least 300 particles within 
L 200 (M 200 ~ 3.4 X 10^ h~^ Mq). This 300-particle limit was 
derived by comparing the mass function of COCO-COLD with that 
of its lower-resolution counterpart COLOR-COLD. Below this limit, 
the results of the simulations become increasingly unreliable. The 
results at high masses are noisy because of the small number of 
high-mass haloes formed in the relatively small volume of our sim¬ 
ulations. 

The general trend across redshifts is similar: for haloes with 
M 200 >2x10® h~^ Mq, the dHMF in COCO-WARM and COCO¬ 
COLD are almost identical. The abundance of haloes below this 
mass scale is strongly suppressed in COCO-WARM, to the extent 
that, at z = 10, there are 5 times fewer ~ 10® h~^ Mq haloes than 
in COCO-COLD. The delayed non-linear structure formation below 
~ 2 X 10® h~^ Mq can also be seen from the fact that there are 
as many haloes with M 200 = 10® h~^ Mq in COCO-WARM at 
2 : = 10, as there are haloes with M 200 = 6 x 10® h~^ Mq in 
CO CO-COLD at that redshift. 



Figure 4. The median redshift of formation of all FOF groups in COCO- 
WARM and COCO-COLD, as a function of the halo mass, M 200 . The red¬ 
shift 2;form is defined in the text. The error bars / shaded region represent 
the bootstrapped errors on the median in each mass bin in COCO-WARM 
and COCO-COLD respectively. As expected, there is good agreement at the 
high-mass end, whereas the differences between CDM and WDM become 
apparent below ~ 2 x 10^ h~^ Mq . The thin red line is a fit to the COCO- 
WARM redshift of formation, using Eq. 8. 


Within the CDM paradigm, there are a number of analytic pre¬ 
dictions for the differential halo mass function (dHMF), notably 
the Press-Schechter formula (Press & Schechter 1974; Bond et al. 
1991; Lacey & Cole 1993), and the ellipsoidal collapse model (ST; 
Sheth & Tormen 1999, although this model is not fully analytic 
since it is tuned to numerical simulations). The dHMF is given by: 

dlogM M'^ ^ ^ 

where f{iy) is the so-called halo multiplicity function and for hi¬ 
erarchical cosmologies has a universal form (see e.g. Jenkins et al. 
2001; Reed et al. 2007; Tinker et al. 2008; Angulo et al. 2012). In 
the ST formalism, it is approximated by: 

/ = ^\/— [1 + (10) 

V 71 

Here, v = 5l{z)/a‘^{M),A = 0.3222,^ = 0.707 and p = 0.3. 
In linear theory, Sc{z) = 1 .686/79(2;), where 79(2;) is the linear 
growth rate of perturbations. The value of Sc is appropriate for the 
Einstein-de Sitter model, but differs slightly in ACDM due to a 
weak dependence on Qrn{z). Finally, a‘^{M) is the variance in the 
mass density field on mass scale, M, given by: 

a^{M) = J ^A^{k)W^{k,M). ( 11 ) 

Here, W{k, M) is the Fourier transform of a window function con¬ 
taining mass M, and A^(/c) is the dimensionless power spectrum 
as defined in Fig. 1. 

In the Press-Schechter and Sheth-Tormen formalisms, the rms 


d log a 


dlogM 


(9) 
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Figure 5. The redshift evolution of the halo mass function in COCO-COLD 
and COCO-WARM. The solid lines show the CDM results, with the shaded 
regions representing the associated la Poisson errors. The dashed lines with 
error bars represent the equivalent relation from COCO-WARM, with spuri¬ 
ous haloes removed. The different colours show results for a selection of 
redshifts, as indicated in the legend. The grey shaded region corresponds to 
haloes with fewer than 300 particles within r2oo • 


fluctuation amplitude, cr^(M), is assumed to be a monotonically 
increasing function of M. This is no longer true for the truncated 
power spectrum of WDM, so care must be taken when choosing 
an appropriate window function. In the CDM, VC(/c, M) is usually 
chosen to be the real-space spherical top-hat function, a choice that 
results in an excellent match to the dHMF in cosmological N -body 
simulations. The same for WDM predicts an excess of low-mass 
haloes compared to simulations (Barkana et al. 2001; Menci et al. 
2012; Schneider et al. 2012, but see also Schneider et al. 2013). 
This problem was solved by Benson et al. (2013), who generalised 
the (extended) Press-Schechter formalism by using the correct solu¬ 
tion for the excursion set barrier first-crossing distribution in WDM 
models. Rather than the top-hat real-space window function, they 
used a sharp /c-space Alter for WDM, so that the variance, cr(M), 
remains flat up to the half-mode mass and then declines with in¬ 
creasing mass (see Fig. 6). In this formalism the smoothing scale, 
R, is defined as: 


where ks — ^nKjoL, a as defined in Eq. 3, = 0.361 and a = 

2.5. Benson et al. (2013) choose the free parameters such that the 
theoretical mass function turns over at the same scale as the halo 
mass function from simulations. This choice of parameters should 
be applicable to all thermal WDM models, since the effect of the 
WDM suppression is captured in the value of a (Eq. 3). 

In Eig. 7, we compare the z = 0 dHME for COCO-COLD 
(blue squares), the full COCO-WARM (genuine and spurious ob¬ 
jects; green diamonds), the spurious COCO-WARM objects only 


Figure 6. The fractional variance of the density field, a{M), calculated in 
Eq. 11 using a top-hat filter in real-space for CDM, and a sharp A:-space 
filter for WDM. The flattening of the relation below 10® h~^ Mq is due 
to the suppression of power below these scales in WDM, relative to CDM. 
The dashed line indicates the upper limit to the halo masses formed in our 
volume-limited simulations. 


(yellow stars) and the genuine COCO-WARM haloes only (red cir¬ 
cles). 

The solid and dashed black lines in Fig. 7 show the ST predic¬ 
tions for the mass functions in CDM and WDM respectively. For 
M 200 > 2 X 10^ h~^ Mq, the mass functions for CDM and WDM 
trace one another exactly, as expected. Below this mass, the WDM 
mass function begins to peel off from the CDM case, reaching half 
the CDM amplitude at M 200 ~ 2 x 10® h~^ Mq. This agrees with 
the half-mode mass scale, Mhm, introduced in Section 2.1. The raw 
WDM mass function (green diamonds) is entirely dominated by the 
spurious objects (yellow stars) below ~ 4 X 10’^ Mq, where 
the mass function shows an artificial upturn. On the other hand, 
the cleaned WDM sample, represented by the red circles, contin¬ 
ues to fall off smoothly from the regime free of artificial haloes. 
The feature at ~ 2 x 10^ h~^ Mq could be related to the cut, 
Mmax = 3.2xl0^/i“^ Mq, applied as part of the cleaning proce¬ 
dure (Section 2.3), but, in any case, this is very close to the resolu¬ 
tion limit which is also the mass scale at which the spurious haloes 
begin to dominate the mass function. 

The main conclusion to be drawn from Fig. 7 is that above 
the resolution limit, the modified ellipsoidal collapse model repro¬ 
duces the WDM mass function remarkably accurately, over nearly 
6 orders of magnitude in mass. 

3.3 Halo density and mass profiles 

Spherically-averaged radial density profiles provide the simplest 
and most direct descriptor of halo structure. We calculate profiles 
in radial shells equally-spaced in log(r/r 2 oo). As we discussed 
in Section 2.2, haloes of mass above 10® h~^ Mq can be bijec- 


,15 


MNRAS 000, 1-16 (2015) 















































Halo properties in the COCO WDM simulation 9 



Figure 7. Differential halo mass functions from the COCO-WARM and 
COCO-COLD simulations, compared to the predictions of the ellipsoidal col¬ 
lapse formalism of Sheth & Tormen (1999). The solid lines show the pre¬ 
dictions of the standard formalism applied to CDM; the dashed lines show 
the predictions of the modified, sharp fc-space filter of Benson et al. (2013). 
The symbols represent results from our simulations as denoted in the leg¬ 
end: blue squares for COCO-COLD, green diamonds for all COCO-WARM 
FOF haloes, red circles for the genuine haloes and yellow stars for spurious 
haloes. The grey shaded region corresponds to haloes with fewer than 300 
particles within r 2 oo- 


lively matched in COCO-WARM and COCO-COLD. To compare den¬ 
sity profiles in the two models, we stack the individual profiles of 
matched and dynamically relaxed haloes in narrow bins of halo 
mass of width A log(M 2 oo) = 0.3. To determine whether or not a 
halo is relaxed, we make use of the criteria for dynamical equilib¬ 
rium set out by Neto et al. (2007): (1) the displacement of the centre 
of mass from the potential centre should be less than O.OTrvir and 
(2) less than 10% of the mass within Lvir should be in the form of 
substructure. 

The stacked differential density profiles are shown in Fig. 8 
for a variety of mass bins, with the ratio of the densities 
shown in the bottom panels. For masses sufficiently larger than 
~ 2 X 10^ h~^ Mq, we expect negligible differences in the 

properties of CDM and WDM haloes: this is apparent in mass bins 
with M 200 > 10^^ h~^ Mq. Systematic differences in the density 
profiles begin to appear at around M 200 ~ 5 x 10^° h~^ Mq: the 
WDM haloes have slightly but systematically lower central densi¬ 
ties than their CDM counterparts. This halo mass is two orders of 
magnitude higher than the half-mode mass, and an order of mag¬ 
nitude higher than the scale at which the mass functions begin to 
differ (Fig. 7). The difference in central density grows as the mass 
decreases and reaches ~ 30% at the smallest mass bin shown, 
M 200 ~ 1-4 X 10^ h~^ Mq. We discuss the physical reason for 
this in the next section. 

It is now well established that the density profiles of dark 
matter haloes in general are well described by the NFW profile 


(Navarro et al. 1996, 1997): 

^ _ (13) 

Pc ir/rs){l+r/rsf ' 

where Sc is a characteristic overdensity and Vs is a scale radius. 
These two parameters are strongly correlated and depend only on 
halo mass (Navarro et al. 1997). The NFW form is a nearly uni¬ 
versal profile in the sense that it approximately fits the profiles 
of relaxed haloes of any mass formed by gravitational instability 
from all the initial conditions and cosmological parameters that 
have been tested so far. The universality of the NFW profile is in¬ 
timately related to the way in which haloes are assembled (Ludlow 
et al. 2013). 

We fit NFW profiles to the stacked density profiles of COCO¬ 
WARM and COCO-COLD in Fig. 8, between the radial range defined 
by the Power et al. (2003) convergence radius, rconv (defined as the 
radius within which the relaxation time is of the order of the age of 
the Universe), and r 2 oo, minimising the following quantity: 

2 A^bins 

^ AT- ZM XI -Inp^Fw {Sc]rs)f . (14) 

z=l 

We obtain the best-fitting values of the scale radius, r^, which de¬ 
fines the halo concentration, C 200 = ^200 /Ls . This parameter pro¬ 
vides a unique characterisation of the NFW density profile; the val¬ 
ues of C 200 for the stacked profiles are quoted in Fig. 8. There is a 
clear trend in that for large halo masses, where the density profiles 
in COCO-WARM and COCO-COLD are similar, the concentrations 
are nearly identical but, for masses below ~ 5 x 10 ^° h-'^ Me, the 
concentrations of WDM haloes are systematically lower than those 
of CDM haloes. 

In many cases, even better fits to the density profile are pro¬ 
vided by a formula first used by Einasto (1965) to describe star 
counts in the Milky Way. This formula, which has an additional 
free parameter, was dubbed the “Einasto profile” by Navarro et al. 
(2004), who showed that it provides a very good fit to CDM haloes: 



where p -2 is the density at r = r- 2 , the radius at which the loga¬ 
rithmic slope of the profile is —2 (or where r^p has its maximum). 
The parameter r _2 in the Einasto profile is analogous to the scale 
radius, r^, of the NFW profile. This allows an equivalent definition 
of halo concentration, C 200 = r 2 oolr- 2 . The parameter a (not to 
be confused with the one in Eq. 3) is a shape parameter that con¬ 
trols the curvature of the profile in the inner regions. A value of 
a 0.17 results in a good match to CDM haloes over a wide 
range of masses (Navarro et al. 2004; Gao et al. 2008). 

This is demonstrated in Eig. 9, which is similar to Eig. 8, 
but with Einasto profiles fitted instead of NEW profiles. It is ap¬ 
parent that the shape parameter, a, allows a better fit to the halo 
density profiles in both COCO-WARM and COCO-COLD, especially 
in the inner parts. It is also interesting to note that the concentra¬ 
tions inferred from the Einasto profile fits tend to be slightly lower 
than those inferred from the NEW profile fits especially at higher 
masses. 

In Eig. 10 we compare the ratio of M 200 values for indi¬ 
vidually matched haloes in COCO-WARM and COCO-COLD at the 
present day. We consider only haloes with M 200 > 10® h~^ Mq 
for which we have almost complete matching (~ 97%) between 
the two simulations, and plot the ratio, as a func¬ 

tion of M §^^. The solid red line shows the median ratio, whereas 
the dashed red lines represent the 16-th and 84-th percentiles. The 
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r/r2oo ^ r/r2oo 




r/r2oo ^ r/r2oo 




Figure 8. Stacked spherically-averaged density profiles in coco-warm (red) and coco-cold (blue). For each mass bin we compare the profiles of only 
relaxed, matched haloes in the two simulations; the number in each bin is indicated in each subpanel. The vertical dashed line represents the convergence 
radius, Vconv, and filled symbols indicate the range of the profile above this limit, whereas open symbols denote the radial range below it. The dashed red and 
blue lines are NFW fits to the WDM and CDM profiles respectively. Note that the density profiles have been scaled by (r/r 2 oo)^ so as to reduce the dynamic 
range on the vertical axis. The bottom panels show the ratio of the WDM and CDM densities in each bin. 


masses are very similar for objects > 5 x 10^° h~^ Mq, where 
the ratios agree to within 1%. For masses lower than this, WDM 
haloes are systematically less massive than their CDM counter¬ 
parts, with the deficit in WDM halo mass reaching ~ 30% at 
= 10^ h~^ Mq. Flaloes of these masses in WDM form 
later than their CDM counterparts and thus have less time to grow. 

In Fig. 11 we show the cumulative radial distribution of 
mass in haloes in COCO-WARM (red lines) and COCO-COLD (blue 


squares). The ratios are shown in the lower panels. From Fig. 10, 
we expect the cumulative profiles to be very similar at r/r 2 oo = 1 
except in the lowest mass bin, where WDM haloes are slightly (~ 
10%) less massive than their CDM matches. The same trend seen in 
the density profiles is apparent here: for M 200 < 5 x 10^° h~^ Mq, 
the profiles are less concentrated in the central regions in COCO¬ 
WARM than in COCO-COLD. The reason for this difference is dis¬ 
cussed in the next section. 
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Figure 9. Same as Fig. 8, but with Einasto fits to the COCO-WARM and COCO-COLD density profiles. 


3.4 The concentration-mass relation 

As mentioned in the previous section, the density profile of a dark 
matter halo is characterised by its concentration. As a result of their 
hierarchical formation process, the inner parts of haloes in CDM 
and WDM are essentially in place even before the bulk of the halo 
mass is assembled (Wang et al. 2011). The concentration reflects 
the mean density of the Universe at the epoch when these inner 
regions are in place and the earlier a halo forms, the higher its con¬ 
centration is (Navarro et al. 1997). 

In Section 3.3, we found that the Einasto profile provides a 
slightly better fit to the density profiles of WDM and CDM haloes 


than does the conventional NEW profile. Eurthermore, Einasto fits 
are less sensitive to the radial fitting range (Gao et al. 2008, but 
see also Ludlow et al. 2013). Eor these reasons, we proceed to de¬ 
rive the concentration-mass relation in our simulations using fits of 
the Einasto profile to the density profiles of individual haloes (not 
the stacks). Again, fitting is performed between the convergence 
radius, rconv, and r 2 oo, while minimising the rms of the fit: 


2 

<^fit — 


-j ^hins 

- - - y 

^bins — 1 

1 = 1 


[InPi - InpEin {p-2]r-2;a)f 


(16) 


To obtain the halo M 200 — C 200 relation we first split the haloes 
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two interesting features of note: firstly, for all redshifts, the down¬ 
turn in the WDM halo concentrations occurs at roughly the same 
halo mass, M 200 ~ 5 x 10 ^° h~^ Mq\ and secondly, at fixed 
mass, the ratio decreases with decreasing redshift. The fact that the 
mass at which WDM halo concentrations begin to peel-off from the 
CDM relation is almost independent of redshift refiects the narrow 
redshift range in which the inner parts of WDM haloes form. 

In COCO-WARM we also find that the evolution of the mass- 
concentration relation over redshift can be approximated using a 
simple functional form motivated by Eq. 8 (see Schneider et al. 
2012 ), with an extra redshift-dependent component: 

WDM / ]\/r \ “^2 

^200 _ ( -| I -^^hm \ / H I \f3(z^ 

x(i + »)'’- or, 

Here, Mhm is the half-mode mass, z is the redshift of interest, 
71 = 60, 72 = 0.17 and /3(z) = 0.026z — 0.04. The predic¬ 
tions of our model are shown in the upper panel of Fig. 12 using 
the thin colour lines. While the model does not fully capture the 
relatively flat relationship at z = 3 and 4 in COCO-WARM, it gen¬ 
erally reproduces the trends in the simulation and provides a good 
fit up to z = 2, over nearly 5 orders of magnitude in halo mass. 


Figure 10. Ratio of halo mass (M 200 ) for all (relaxed and unrelaxed) 
matched haloes above M 200 > 10® h~^ Mq in coco-warm and coco¬ 
cold, as function of The solid red line shows the median relation 

in bins of whereas the dashed red lines indicate the 16-th and 84th 

percentiles. 


into bins equally spaced in logarithmic mass. We then fit an Einasto 
profile to each halo individually, removing all unrelaxed haloes ac¬ 
cording to the Neto et al. (2007) criteria. We then find the median 
halo concentration in each mass bin and estimate its uncertainty 
using bootstrap resampling. 

Fig. 12 shows the (median) concentration-mass relations for 
COCO-COLD (dotted lines and shaded regions) and COCO-WARM 
(points with error bars) at redshifts z = 0,0.5,1,2,3 and 4 (differ¬ 
ent colours as indicated in the legend). These relations display the 
same qualitative behaviour seen in the density profiles in Fig. 11. 
For haloes with mass M 200 > 10^^ h~^ Mq, the concentra¬ 

tions of CDM and WDM haloes agree well over all redshifts. For 
masses below this value, WDM haloes are less concentrated than 
their CDM counterparts at all redshifts. This is a direct result of 
the later formation epoch of haloes of a given mass in WDM, and 
refiects the fact that the mass within r _2 in WDM haloes is assem¬ 
bled when the background density of the Universe is lower than in 
the CDM case. 

Whereas the CDM halo concentrations continue to increase as 
power laws towards lower masses, reflecting hierarchical growth, 
the WDM halo concentrations turn over at M 200 < 5 x 

10^° h~^ Mq and eventually begin to decrease (see also Schneider 
et al. 2012; Maccio et al. 2013). This echoes the finding in Fig. 11 
that the mass in the central regions of WDM haloes begins to fall 
below that in the CDM case roughly below this mass. This mass 
is an order of magnitude larger than the mass scale at which the 
mass functions begin to differ (~ 2 X 10® h-^ Mq, see Figs. 4, 7). 
This result is not entirely surprising: the concentration is sensitive 
to the inner parts of the profile and it is this inner mass (which we 
can roughly identify with the matter contained within r- 2 ) which 
is assembled later in WDM than in CDM, while most of the mass 
actually lies in the outer parts of the halo. 

The lower panel of Fig. 12 shows the ratio of the concentra¬ 
tions in COCO-WARM and COCO-COLD, c^™/c§)™. There are 


3.5 The shapes and spins of haloes 


In this section we examine the shapes and spins of WDM haloes. 
The shapes are most commonly quantified by the triaxiality, de¬ 
fined through the halo inertia tensor: 

N 200 

lij — TTlp ^ ^ Xn,iXn,j 5 ( 1 ^) 

n = l 


where A^ 2 oo is the number of particles within r 2 oo, rrip is the mass 
of the simulation particle, and Xn,i is the ith coordinate of the 
n—th particle relative to the halo centre. The eigenvalues of the 
inertia tensor define the axial lengths of an equivalent uniform den¬ 
sity ellipsoid, a ^ b ^ c, which can be related to those of the 
halo itself (Bett et al. 2007). The sphericity is defined as c/a (as in 
Section 2.3): the higher its value, the less aspherical the ellipsoid’s 
projection. The triaxiality is defined as T = (a^ ~ ~ c^): 

large values correspond to prolate ellipsoids, small values to oblate 
ellipsoids. 

The results for our simulations are shown in Fig. 13, where 
blue represents CDM and red WDM, with the top panel compar¬ 
ing the median triaxiality, and the lower the median sphericity. Er¬ 
rors on the median quantities were obtained by bootstrap resam¬ 
pling. Previous iV-body simulations of CDM haloes have shown 
that triaxiality correlates with halo mass, with triaxility decreasing 
with decreasing halo mass (Frenk et al. 1988; Allgood et al. 2006; 
Munoz-Cuartas et al. 2011; Maccio et al. 2013). This trend refiects, 
in part, the younger dynamical age of more massive haloes. Fig. 13 
shows that the same trend is present for WDM haloes but below 
M 200 ~ 10^° h~^ Mq, WDM haloes are slightly less triaxial than 
their CDM counterparts. 

A more significant trend is revealed when comparing the spin 
of haloes in the two simulations. The spin is best characterised by 
the parameter. A, defined as: 


jywi 

GM^/^ 


(19) 


(Peebles 1969), where J is the magnitude of the angular momen¬ 
tum of the halo, E is its total energy and M is the mass (which 
we take to be M 200 ). Haloes acquire a net angular momentum 
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Figure 11. Stacked cumulative mass profiles of relaxed, matched haloes in different mass bins for WDM (solid red lines) and CDM (blue squares). The 
lower panels show the ratio of the WDM mass to the CDM mass as a function of radius from the centre of the halo (in units of r 2 oo)- For haloes with 
^200 > 10^^ Mq, the mass profiles are nearly identical, but below M 200 ^ 5 x 10^° h~^ Mq they differ noticeably. 


through tidal torques during growth in the linear regime which 
can be subsequently modified and rearranged by mergers (Peebles 
1969; Doroshkevich 1970; White 1984). Since the merger histories 
are different for CDM and WDM haloes, we might expect some 
differences in their final angular momentum configurations. In par¬ 
ticular, given that tidal forces associated with mergers tend to re¬ 
distribute angular momentum from the central parts of haloes to the 
rest of the halo, the smaller frequency of mergers in WDM might 


facilitate the formation of extended spinning galactic disks (Frenk 
et al. 1988; Navarro & Benz 1991; Navarro & White 1994) 

The spin parameters in our two simulations are compared in 
the top panel of Fig. 14. Previous cosmological CDM simulations 


^ We note that the inability of many early simulations to form extended 
disks in the CDM model - the so-called “angular momentum” problem - 
is readily solved when appropriate prescriptions for supernovae feedback 
are included in the simulations (see e.g. Okamoto et al. 2005; Scannapieco 
etal. 2011). 
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Figure 12. The median concentration-mass relation and its redshift evolu¬ 
tion for haloes in COCO-COLD and COCO-WARM. The colour dotted lines 
show the median relation over redshift for CDM haloes, as indicated in the 
legend. The shaded regions represent the errors in the median, as estimated 
by bootstrap resampling. The points with the error bars show the corre¬ 
sponding redshift relation in WDM. Only relaxed haloes are included. The 
thin colour lines show the results of the fitting formula introduced in Eq. 17. 


showed a very weak correlation between spin and halo mass, with 
a median value of A 0.033, across a wide range of halo masses 
(Davis et al. 1985; Barnes & Efstathiou 1987; Warren et al. 1992; 
Steinmetz & Bartelmann 1995; Cole & Lacey 1996a,b; Mo et al. 
1998; Bett et al. 2007). Our COCO-COLD simulation reproduces 
this trend and extends it to lower masses, M 200 = 10® h~^ Mq. 

For M 200 > 5 X 10^° h~^ M©, the A values for WDM haloes 
are almost identical to those of their CDM counterparts. However, 
for smaller halo masses A decreases systematically with decreas¬ 
ing mass and is lower than the CDM value by almost 30% at 
M 200 ~ 10® h~^ Mq. This is consistent with the results of Bul¬ 
lock et al. (2002), who found that three out of four haloes below 
the WDM cutoff in their simulation had lower values of A than the 
equivalent CDM matches. Note that in the top panel of Fig. 14 we 
include all haloes, not necessarily matches, which explains why in 
some of the largest mass bins, the median spins are not exactly the 
same in WDM and CDM. In addition, we only include haloes with 
more than 1000 particles within r 2 oo since particle shot noise dom¬ 
inates the estimates of angular momentum for low particle numbers 
(Frenk et al. 1988; Bett et al. 2007, although we use a more conser¬ 
vative lower limit than the latter’s choice of 300 particles). 

To investigate why the spins of dwarf galaxy haloes are lower 
in WDM than in CDM we consider the relative contributions of 
energy, angular momentum and M 200 to A, illustrated in the bot¬ 
tom panel of Fig. 14, this time for bijectively matched haloes. 
The ratio of the median spin parameters is shown by the black 
squares and the ratio of the geometric means of the quantities that 
enter into Fq. 19 are shown by the other colour lines (magenta 
for Jcdm/Jwdm, cyan for |F;cdm/L^wdm|^'^^ and yellow for 



Figure 13. Median halo triaxiality (top panel) and halo sphericity (lower 
panel) in COCO-WARM (red points) and COCO-COLD (blue lines). The er¬ 
rors on the median were obtained by bootstrapping 100 different samples in 
each case, and is represented by the red error bars for WDM and the blue 
shaded region for CDM. Only particles within r 2 oo were used to compute 
these properties from the inertia tensor. 

(M 2 oo,cdm/M 2 oo,wdm)^^^). The combination of these ratios in 
Fq. 19 should reproduce the ratio of spin parameters, and this is 
shown in the thick green line. Part of the reason for lower WDM 
spins below ~ 10^° h~^ Mq is their slightly lower total energy 
which results from their lower concentration. The dominant fac¬ 
tor, however, is their lower angular momentum relative to CDM 
haloes, ~ 25% at 10® h~^ Mq. The cause of this could be related 
to the differing merger histories in WDM and CDM and the likely 
more quiescent mass accretion of WDM haloes which can result in 
smaller spins (Bullock et al. 2002; Vitvitska et al. 2002; Hetznecker 
& Burkert 2006). 


4 CONCLUSIONS 

We have presented results from the COpernicus COmplexio 
project, a set of cosmological “zoom” simulations in which the dark 
matter is assumed to be either CDM (COCO-COLD) or a thermal 3.3 
keV WDM particle (COCO-WARM). The combination of mass res¬ 
olution and volume of our simulations provides a rich statistical 
sample of haloes over seven decades in mass. This WDM model 
is particularly interesting because it corresponds to the “warmest” 
particle allowed by current Lyman—a constraints (Viel et al. 2013) 
and has a linear power spectrum cutoff similar to that for the 
“coldest” 7 keV sterile neutrino, evidence for which has recently 
been claimed to be found in galaxies and clusters (Bulbul et al. 
2014a; Boyarsky et al. 2014a). This cutoff - manifest in haloes 
of M 200 ^ 2 X 10^ h~^ Mq for our assumed particle mass - is 
reflected both in the population statistics and the structure of indi¬ 
vidual haloes. 
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Figure 14. Top panel: the median halo spin-mass relation at 2 ; = 0 for 
COCO-WARM (red points) and COCO-COLD (blue line). Errors on the me¬ 
dian for the WDM (shown by error bars) and for CDM (shown as the shaded 
region) haloes were calculated by bootstrap resampling. Bottom panel: the 
relative contributions of energy, angular momentum and halo mass to the 
spin of the halo. The black squares show the ratio (CDM to WDM) of the 
median spin parameters (from the top panel). The magenta, cyan and yellow 
lines show the ratios of the geometric means of the angular momentum, en¬ 
ergy and M 200 respectively, which when multiplied together appropriately 
yield the thick green line, which show the ratio of the geometric means 
of '^CDM and Awdm- As expected, the squares trace out the ratio of the 
geometric means. Note that ratios of all quantities are taken between the 
bijectively matched COCO-WARM and coco-cold haloes. 

The formation of structure begins significantly later in COCO¬ 
WARM than in COCO-COLD. Across all redshifts, differences in the 
halo mass function between COCO-WARM and COCO-COLD begin 
to appear at a mass roughly one order magnitude larger than the 
nominal half-mode mass. Below ~ 2 x 10® h~^ Mq, the WDM 
mass function declines rapidly but there are still some small haloes 
present at suprisingly large redshifts: at z = 10, for example, there 
are almost 5 times as many haloes with M 200 ~ 10® h~^ M© in 
COCO-COLD than in COCO-WARM. We find that the z = 0 halo 
mass functions in both COCO-WARM and COCO-COLD are well 
described by previous analytic fits to the CDM halo mass func¬ 
tion (e.g. Sheth & Tormen 1999) down to our resolution limit, 
M 200 ~ 3 X 10^ h~^ Mq, provided that the window function used 
to compute the mass variance, in the WDM case is cal¬ 

culated using a sharp /c-space filter, as described by Benson et al. 
(2013). 

Just as for COCO-COLD, the spherically averaged density pro¬ 
files of haloes in COCO-WARM, down to dwarf galaxy scales, are 
well described by NEW or Einasto profiles. The concentration- 
mass relation, M200 — C200 (where we have defined concentra¬ 
tion using the Einasto profile), in COCO-WARM begins to peel 
off from the corresponding relation in COCO-COLD at a mass of 
~5xl0^®/i“^ M©, refiecting the later formation epoch of haloes 
of a given mass in WDM compared to CDM. This mass is larger 


than the scale below which the WDM mass function is suppressed 
because halo concentration is determined by the epoch when the 
inner regions of a halo form. The mass at which the concentration 
begins to differ in the two simulations is almost independent of red- 
shift out to z ~ 4. At the present day, the typical concentration of 
a halo of mass 10® h~^ Mq in COCO-WARM is C200 — 8 com¬ 
pared to C200 — 12.7 in COCO-COLD. The trends and evolution of 
the concentration-mass relation can be approximated by the fitting 
formula provided in Eq. 17. 

The generally triaxial shapes of haloes in COCO-WARM and 
COCO-COLD are very similar. However, we find that, for masses 
below ~ 5 X 10^® h~^ Mq, WDM haloes have slightly lower val¬ 
ues of the spin parameter. A, (up to 30%) than their CDM coun¬ 
terparts. Whether this and the other differences in halo properties 
between WDM and CDM that we have found have implications 
for galaxy formation needs to be investigated using hydrodynamic 
simulations. 
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